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We studied the energy levels of graphene based Andreev billiards consisting of a superconductor 
region on top of a monolayer graphene sheet. For the case of Andreev retro-reflection we show 
that the graphene based Andreev billiard can be mapped to the normal metal-superconducting 
billiards with the same geometry. We also derived a semiclassical quantization rule in graphene 
based Andreev billiards. The exact and the semiclassically obtained spectrum agree very well both 
for the case of Andreev retro-reflection and specular Andreev reflection. 
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In the well-known Andreev billiards consisting of a 
normal metal surrounded by a superconductor (NS) the 
dynamics of the quasiparticles is determined by the so- 
called Andreev retro-reflection []]]. The spectrum of An- 
dreev billiards is described by the Bogoliubov-de Gennes 
(BdG) equation and has been long studied 0, 0] (for re- 
view of the topic see, eg, [H). 

The electronic properties of graphene can be de- 
scribed accurately by massless Dirac fermion type ex- 
citations using two dimensional relativistic quantum me- 
chanics 0, 0] and also by semiclassical methods [§] 
(for reviews on the physics of graphene see, e.g., [§]). In 
the seminal paper by Beenakker [10| it has been shown 
that when monolayer graphene is interfaced with a su- 
perconductor then two types of Andreev reflection are 

(G) 

possible depending on the ratio of the Fermi energy Ep 

and the electron energy E. For E F G ^ > > E the 

Andreev retro-reflection is dominant as in NS billiards 
(here A( G ) is the superconducting pair potential induced 

in the graphene). In contrast, when Ep < E < A^- G \ 
a different type of scattering process takes place at the 
graphene-superconductor interface, which is named spec- 
ular Andreev reflection . The specular Andreev reflection 
does not exist in NS systems and it is a prominent con- 
sequence of the peculiar band structure of the monolayer 
graphene. Beenakker's paper has been followed by nu- 
merous works (for a review on Andreev reflection in 
graphene see article [HI)- Note that although graphene 
itself is not superconducting, due to the proximity effect 
a superconductor can induce non-zero pair-potential in 
the graphene as well. Indeed, supercurrent has been ob- 
served experimentally [l3| between two superconducting 
electrodes on top of a graphene monolayer. Moreover, ex- 
perimental results of Ref. 1J] attest to the ballistic prop- 
agation of quasiparticles in graphene-superconductor hy- 
brid structures. 



The most widely studied theoretical model of Andreev be linearized around the Fermi energy E^> as E(p) = 



billiards is that of a two dimensional electron gas (2DEG) 
in a quantum dot contacted by a bulk superconductor 
(see eg. [4j). One of the major obstacles that has 
thwarted so far the direct comparison of the theoretical 
predictions and experimental results is the inevitably ex- 
isting tunnel barrier and mismatch of the Fermi- velocities 
and effective masses between the 2DEG and the super- 
conductor (often referred to as "non-ideal NS interface" 
in the literature). This mismatch causes the probability 
of normal reflection to increase at the NS interface while 
the probability of the Andreev reflection diminishes sig- 
nificantly. The situation when both normal and Andreev 
reflection take place at the NS interface is theoretically 
more difficult to address. In graphene however, when the 
superconductivity is induced by external superconduct- 
ing contacts, such mismatch may not exist so that the 
graphene-superconductor systems may experimentally be 
ideal to study most of the theoretical predictions made 
assuming perfect (ie with no mismatch) NS-interfaces. 

In this paper we consider graphene Andreev billiards 
(GABs). In particular, we assume that in a closed region 
T> of the graphene sheet the superconducting pair poten- 
tial is zero and outside this region it takes on a constant 
value A( G ). We demonstrate, in one hand, that when the 
retro-reflection is the dominant scattering process at the 
normal graphene-superconductor interface the electronic 
properties of GABs can indeed be obtained in semiclassi- 
cal approximation from the known results for NS billiards 
with ideal NS interface. On the other hand, we also cal- 
culate the exact spectrum of a GAB for the case when 
the dominant scattering process is the specular Andreev 
reflection and we show that it can also be understood 
using semiclassical considerations. 

To see the relation between the energy spectrum of NS 
billiards and GABs note the following: the dispersion 
relation of the quasiparticles in the normal (A^ = 0) 
region of the NS billiards for energies E < AW can 
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where the sign 



refers to the 



electron-like (hole-like) quasiparticles. Here p is the mag- 



nitude of the momentum, p F 



(N) 



2mE F ^ is the Fermi 

momentum and v F " ' = Pp^/m is the the Fermi veloc- 
ity. This linearization is allowed if we are interested in 
the properties of the bound states (E < of NS 

billiards because for typical NS billiards the dimension- 
less parameter A^) / E^ <C 1 is much less than unity. 
The same linear dispersion can be found for electron- 
like (hole- like) quasiparticles in the A^ G ) = region 
for GABs in the retro-reflection regime but with Fermi 

(G) • (G) 

velocity v F and Fermi momentum p F . This simple 
observation is the core of the intimate relation between 
the graphene based and normal metal Andreev billiards. 
As long as the effect of the superconductor in semi- 
classical approximation can be described by the same 
way for GABs as for the NS billiards, i.e. by a sim- 
ple phase shift — arccos^/A^' )), one can expect that 
when the Andreev retro-reflection is the dominant scat- 
tering process the gross features of the energy spectrum 
of a GAB will closely resemble the spectrum of a NS bil- 
liard having the same geometry. This happens because 
the quasiparticles have linear dispersion in both cases. 
Moreover, if the Fermi velocities and Fermi momentums 



are the same i.e., vp 



hk F ^/m 



E F G) /(hk^') and 



PF 



= \ 2mE 
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E F /vf the quasiparticles in the 
^(n.g) _ q re gj on w \\\ have the same dispersion rela- 
tion for both NS billiards and GABs. Note that if pp 
and v-p are the same then E F N ^ = E F G ^ /2. 

To demonstrate the idea discussed above we consider 
a simple, circular shape GAB. It consists of normal 
graphene region of radius R surrounded by supercon- 
ducting graphene. Owing to the valley degeneracy of 
the Hamiltonian the full BdG equation for graphene- 
superconductor systems decouples to two four by four, 
reduced Hamiltonians that are related to each other by a 
unitary transformation (see, e.g., 0). We now take the 
one corresponding to the valley K. Due to the circular 
symmetry of the setup the reduced Hamiltonian is sep- 
arable in polar coordinates and therefore the eigenfunc- 
tions can be labelled by an integer number m correspond- 
ing to the angular momentum quantum number. One can 
show that the ansatz for the wave functions satisfying the 
Schrodinger equation for the reduced Hamiltonian in the 
region where A' G ) = 0, ie, for r < R with energy E 

are 9^{r,<p) = (cf X ( + W) + ^ X ( -\r,<fi))e im ^ 
where x+\r,ip) = [-ij m (k+r), J m+1 (k+r)e iv , 0, 0] T 



and X L N) (r,^) = [0, 0, -%J m {k-r), J m +i(fc_r)e i *'] T are 

, F ±E)/(hv F ). 



(e f ^±e) 



the two eigenstates, and k± 

In the superconducting region r > R where the 
pair potential is A^ G ^ the wave function has the 

form ¥i ] (r, ip) = [cf xf (r, <p) + cL S) X- S) (r, ¥>)) e im *, 



where x+ ( r > 9 s ) 



(S) (S) 



,(S) 
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,( s )/ 



-iH^{q+r),H^ +1 {q + ry* 



The eigenstate 



X- ( r i f) is obtained by the replacement + — * — and 
the first kind of Hankel functions to the second one 



and q± = [e f C) ± i\/[A( G )] 2 - £ 2 J /(hi*), while 7± = 

e ±larccos ( £/A(G) ). Here J m (x) and H*fc 2) {x) are the 
Bessel and the Hankel functions [HI]. To ensure that 
the wave function of the bound states is normalizable, 
the wave function in the superconducting region must go 
to zero as r — ► oo. This condition can be satisfied by 
choosing the appropriate Hankel function in the eigen- 
states x± ( r ) f) [HI • Finally, the unknown coefficients 



and c'f 1 can be determined from the boundary con- 



„(S) 



ditions ^m\r 



(s) 

R,<p) — \&m \T — R,f) valid for any 
ip. Thus, the condition for non-trivial solutions of the 
coefficients and cf can be found from the zeros of a 
four by four determinant. After some algebra we obtain 
a quite simple secular equation for the energy levels with 
fixed angular momentum index m: 

Im { 7+ Dg(m, E) D^(m, E)} = 0, (la) 



H^ ] (q + R) 



Jm(k+R) 

J m+1 (k + R) H^Uq+R) 



(lb) 



and D { GS ){m,E) 



D, 



(+)/ 
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rn 

0,±1,±2 



H£>(q+R) 



GS (m, —E) , and Im{.} and 
* stand for the imaginary part and the complex 
conjugation, respectively. Note that Hm'(q-R) = 
The solutions of Eq. Jl} for m = 
are the exact energy levels of a circular 
shape GAB. Note that Eq. (QJ is valid both in the case 
of Andreev retro-reflection (E F G ^ > A' G ) > E) and for 
specular Andreev reflection (E F G ^ < E < A^). One 
can also notice that the eigenenergies depend only on 
two dimensionless parameters: E F QS> /A^ and ^i G ^ / R, 
where t;^ = hv F G ^ / A( G ' is the coherence length in the 
superconducting graphene. 

We now compare the density of states (DOS) g(E) = 
Snm — E nm ) of a circular shape GAB and of the 
corresponding NS billiard. For details of the calculation 
see a similar calculation for NS billiards in Ref . [l6[ ■ It is 
more convenient to plot the integrated DOS, namely the 
so-called step function N(E) = J2 nm Q{E — E nm ), where 
&(x) is the Heaviside function. Our numerical results for 
j^W/AW = E ( F G) /{2A^) are shown in Fig.Q] One can 
see that the step functions for the considered NS billiard 
and GAB are indeed very similar. It is also clear from 
Fig.Q] that the °OS g(E) = dN(E)/dE shows singulari- 
ties at certain energies En . Singularities of this kind 
arise in the case of NS billiards as well (see e.g. Ref. fl6j ) 
and we shall discuss their origin below. 
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FIG. 1: (color online) The exact step function N(E) for cir- 
cular shape GAB (red line) and NS (blue line) billiards in the 
case of Andreev retro-reflection. The parameters for GAB 
and NS billiard are ^ G) /A (G) = 10 and £^ N) /A (JV) = 5, 
respectively, and £c G ' /R = 0.12 for both cases. The insets 
show the enlarged parts of the main frame. The trivial factor 
2 owing to the valley degeneracy in graphene is not included. 



We now demonstrate on the example of circular GABs 
that at semiclassical level the results for NS billiards and 
GABs can be mapped into each other by choosing the pa- 
rameters appropriately. In numerous works [H,[l7, lf| 1|| 
UIHl it was shown that for NS billiards in semiclassical 
approximation the step function reads 



N BS (E) = M I 1 ~ FK(E)]} > (2a) 

n=0 

, , mr + arccos (£7AW) tm , 

*n(E) = — 4^ L & N) . (2b) 



Here M is the number of open channels in the nor- 
mal region, £c — hv 1 ^^ j A^ is the coherence length 
in the NS system, F(s) = J Q P(s') ds 1 is the inte- 
grated path length distribution and P(s) is the clas- 
sical probability that an electron entering the billiard 
at the NS interface returns to the interface after a 
path of length s. The path length distribution P(s) 
is normalized to one, i.e., P(s)ds = 1 and one 
can see that it is a purely geometry-dependent func- 
tion. In particular, for circular billiards it was found 
that P(s) = ^ 7 =JL—e(2R - s) and M = 

2Trk { F N) R [16]. Finally, the quantity s n (E) in Eq. (|2b|) 
depends on the quantization condition for the periodic 
motion of the electron-hole quasiparicles [13, H3|- As 
it has been pointed out in the introduction, in good 
approximation the quasiparticles have linear dispersion 
in the non-superconducting region for both GABs and 
NS billiards. If the effect of the superconductor in 
GABs can be taken into account by a simple phase shift 
— arccos (E/A^) 0, expressions of the type of Eq. ([2|) 
can be used to calculate the semiclassical approximation 
of N(E) for GABs as well. 



Moreover, employing the same steps as in Ref. [16j . 
from Eq. {1} one can derive the following semiclassical 
quantization rule for circular shape GABs: 



S+ (E) - n r S_ (E) - 2 arccos 



2tt hi 



S±(E)=2d(\k±\R) 



E 



2 I ml arccos 



k±\R' 



(3a) 



(3b) 



where fi r = 1,-1 for Andreev retro-reflection and spec- 
ular Andreev reflection, respectively, and n is a non- 
negative integer. Functions S±(E) are the radial action 
(in units of H) of electrons and holes [22j and the term 
—2 arccos E/A^ in Eq. (f3"aj) accounts for the two An- 
dreev reflections in one period of the orbit, while the 
second term in the left hand side of Eq. (|3"a|) results from 
the sum and the difference of the Maslov indices ir/4 of 
the electron-like and hole-like particles for fi r = 1 and 
fi r = — 1, respectively. 

Formally, in the case of Andreev retro-reflection the 
quantization condition shown in Eq. §3§ is the same as 
for a circular NS billiard [l6j |. but the meaning of k± 
is different for the two systems (for NS billiards see eg. 
Ref. [Hj]). However, from Eq. ^ it is easy to find that if 
HW/fW = p(G)/^f ) and 4 Ar )/A(JV) = 4 G )/(2A(G)) 
then to first order in E/ A^ N ' G ^ the quantization condi- 
tion for circular GABs and NS billiards is the same and 
the step function N(E) is given by Eq. ([2]) with coher- 
ence length ^c G ^ ■ The exact and semiclassically calcu- 
lated N(E) are plotted in Fig.O The agreement between 
the two results is excellent. Moreover, from Eq. @, we 
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FIG. 2: (color online) The exact (red line) and the semiclassi- 
cally (blue line) calculated step function N(E) obtained from 
Eqs. (TfJ) and ([2j|, respectively for the case of Andreev retro- 
reflection. The parameters are the same as in Fig. ([TJ. The 
insets show the enlarged parts of the main frame. 

find that the positions of the singularities in the DOS 
are given by 4 si " s) /A (G) = (n + 1/2) tt/(1 + 2i?/^ G) ) 
valid for such integers n that 4 si " s) < A( G ' holds. Note 
that the position E { n' m&) / A^ of the singularities de- 
pends only on R/^i G ^ but not on 

4 G) / A(G) - Therefore 
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even if E^/A^ ? Ef ) /{2A^) but RW/£W = 

R^/d G \ the singularities in the DOS for a circular 
GAB and NS billiard would appear at the same energies. 

Next we consider the case of specular Andreev reflec- 
tion in graphene Andreev billiards. Again, the solutions 
of Eqs. ([1]) and ([3]) give the exact and the semiclassi- 
cally calculated energy levels of circular shape GABs. In 
Fig-Elthe calculated step function N(E) is plotted and as 
one can see it is completely different from that obtained 
for the case of Andreev retro-reflection shown in Fig. [T] 
Moreover, the quantum results in Fig. [3] again show very 
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FIG. 3: (color online) The exact (red line) and the semiclassi- 
cally (blue line) calculated step function N(E) obtained from 
IU) and ((3]), respectively for specular Andreev reflection. The 
parameters are E^ = and ^i G ^ /R = 0.03. The insets show 
the enlarged parts of the main frame. 

good agreement with the semiclassical ones that can be 
obtained from Eq. ([3|). This implies that in the case of 
specular Andreev reflection the DOS depends linearly on 
the energy for E — > 0. Namely, it can be shown from 
Eq. that in this limit the DOS in semiclassical ap- 
proximation (without the valley degeneracy) is given by 
p(E) = S ^^^yi , where A is the area of the billiard. It 
is interesting to note therefore that (apart from the val- 
ley degeneracy) p{E) is bigger by a factor of 16/-7T 2 than 
in the case of neutrino billiards 23j. 

In summary, we calculated the energy levels of 
graphene based Andreev billiards. We showed that for 
energy levels corresponding to the case of Andreev retro- 
reflection the graphene based Andreev billiards in a very 
good approximation can be mapped to the normal metal- 
superconducting billiards with the same geometry. We 
also derived a semiclassical quantization rule in graphene 
based Andreev billiards and the spectrum obtained from 
this rule agrees very well with that obtained from the 
exact quantum calculations for circular shape of GS bil- 
liards. 
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